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A Maxwell Construction for Phase Separation in Vibrated Granular Matter? 
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Experiments and computer simulations are carried out to investigate ordering principles in a gran¬ 
ular gas which phase separates under vibration. The densities of the dilute and the dense phase are 
found to follow a lever rule. A Maxwell construction is found to predict both the coexisting pressure 
and binodal densities remarkably well, despite the fact that the pressure-volume characteristic P(v) 
is not an isotherm. Although the system is far from equilibrium and energy conservation is strongly 
violated, we derive the construction from the minimization of mechanical work and fluctuating 
particle currents. 


Many-particle systems driven far from equilibrium, 
which occur abundantly in nature, technology, as well as 
in laboratory settings, often exhibit remarkable collec¬ 
tive behaviour m- In spite of the importance of such 
phenomena, the search for underlying principles govern¬ 
ing their dynamics and emerging patterns is still con¬ 
tinuing [Mg. Inspired by analogous problems in equi¬ 
librium thermodynamics, it has proven useful to study 
non-equilibrium steady states (NESS) which are charac¬ 
terized by time-independent macroscopic quantities and 
their fluctuations. 

A paradigmatic model system exhibiting a NESS is 
a driven granular gas [I3fll6| . In its simplest form, a 
granular gas is a cloud of dissipative spherical particles 
maintained in a steady state by a continuous external 
drive HZ1 HU. Recent work has demonstrated that con¬ 
fined grains driven by a periodic external force can sep¬ 
arate into a dilute phase and a dense phase via spin- 
odal decomposition mm, and exhibit behavior similar 
to the phase separation in a van der Waals gas I2TI[121| . 
In an equilibrium fluid the pressure at coexistence, P*, 
as well as the binodal densities (j) g (gas) and <f>i (liquid) 
are determined by the requirements of mechanical and 
thermal equilibrium and by the minimization of the ap¬ 
propriate free energy. These constraints give rise to the 
Maxwell construction on the non-monotonic pressure- 
volume isotherm P(v) [25] . This equal-areas construction 
uniquely identifies P* = Pea , the equal-areas pressure. 

In a granular gas, the degree of dissipation is quanti¬ 
fied by the restitution coeffitient, e < 1, which denotes 
the ration of the relative momenta after/before impact of 
two particles. Clustering is observed whenever e < 1. If 
e is close to unity (elastic limit), it is not surprising that 
the overall clustering behaviour can be described by con¬ 
cepts borrowed from equilibrium statistical physics [23- 
[251125] . However, as the system is taken far away from 
equilibrium, its behaviour is expected to differ qualita¬ 
tively from its equilibrium counterparts. In particular, 
a Maxwell construction, which directly derives from the 
minimization of a free energy functional, is not expected 
to hold. Surprisingly, we find that an equal-areas rule 
(Maxwell construction) quite accurately predicts the co¬ 
existing binodal densities and the pressure for the liquid- 


gas phase separation to a few percent, even if the system 
is well remote from the elastic limit (with e down to 0.65). 
We demonstrate that this results from the minimization 
of residual mechanical work associated with the fluctua¬ 
tions in the system. 

Our experimental apparatus is very similar to that 
used previously [25] , Glass particles were confined be¬ 
tween parallel horizontal plates in a long, thin cell. The 
particles were sieved and selected under a microscope to 
obtain a sample of approximately monodisperse spheres 
with diameter d = 610 /mi. The cell was constructed from 
a lower plate of 3 mm thick, anodized aluminum and a 
top plate of 3 mm thick glass. The plates were separated 
by 10 mm high aluminum walls which also confined the 
particles horizontally such that the internal length, width 
and height of the cell were 280 mm, 10 mm, and 10 mm, 
respectively. The cell was driven sinusoidally in the ver¬ 
tical direction, with variable amplitude, A. The driving 
frequency was kept fixed at 60 Hz. The mean density is 
defined as </> = Nv p /V, where N is the number of par¬ 
ticles, v p is the volume of a single particle and V is the 
internal volume of the cell. Care was taken to ensure that 
the cell was level prior to each experimental run. 

Density profiles measured in experiments using a long 
cell are shown in the main panel of Fig. [T] the corre¬ 
sponding interface positions are shown in the inset. As 
(f> increases, the volume of the liquid phase increases, 
moving the interface to the left. The densities (f>i and 
< j) g appear to be independent of cf). In the inset the lin¬ 
ear fit demonstrates that the system obeys a lever rule, 
(j) = ipifi + <j)gf g , where fi is the liquid fraction and f g 
is the gas fraction. The lever-rule confirms that there 
is an intrinsic mechanism which selects the liquid and 
gas densities as intensive quantities of the system. In 
thermodynamic equilibrium, the steady state is defined 
by maximizing the entropy of the system subject to the 
conservation of energy and momentum. In the NESS, 
energy conservation is violated due to the external driv¬ 
ing and the dissipation of energy into internal degrees of 
freedom of the particles. This invites the question: what 
sets the coexisting densities in this NESS? 

In addition to our experiments, we have also carried 
out time-driven molecular dynamics simulations. The 
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FIG. 1. Main panel: the mean grey level from photographs 
of the experiment. Data are shown for densities in the range 
0.04 < (j> < 0.12, increasing from right to left. The driv¬ 
ing amplitude is A=2. Id. The local density, </>(*), is inferred 
from the average gray level in top-view photographs of the 
system, and scaled such that the peak values match the val¬ 
ues predicted by the lever rule plot by extrapolation to fi=0 
and fi = 1. Inset: The crosses mark the interface position ob¬ 
tained from the experiment. The dashed line is linear fits to 
the data. Zero liquid fraction corresponds to gas density, 4> g , 
unity liquid fraction represents liquid density, (pi (circles). 



FIG. 2. Dependence of the pressure, P, on the volume per 
particle, v. The solid line shows the pressure in a small 
square-base cell with side length L=20d. In the unstable re¬ 
gion, the pressure tends towards the dashed tie line which 
connects the binodal densities calculated in the large cell 
(Z/=460d). The open circles, squares and diamonds show 
the pressure calculated in cells of length L=80d, 120d, 160d, 
respectively, demonstrating the convergence to the large-cell 
limit. The hatched areas illustrate the equal-areas rule. The 
grey, dashed asymptotes schematically indicate the low tem¬ 
perature branch (left) and high temperature branch (right). 


simulations have previously been shown to accurately 
capture the physics of the system under study [20;, |24] . 
The particles are modelled as monodisperse soft-spheres 
with diameter d = 610 /im. Dissipation is included by 
a normal coefficient of restitution, e (implemented using 
a linear-spring and dash-pot damping) [27]. The effects 
of tangential forces and rotational degrees of freedom are 
neglected because these have been shown to have mi¬ 
nor impact on the physics of the system [2HH3D3 - The 
simulated cell has dimensions 460d x 20d x 16. Ad, thus 
closely resembling the experimental system. Reflecting 
boundary conditions on the short walls of the cell are 
used in order to study a single interface between the two 
coexisting phases. On the long walls, periodic boundary 
conditions are employed. The results presented in this 
paper are for e = 0.8, matching well with the experi¬ 
ment, but similar results are found for the rather wide 
range of 0.65 ^ e ^ 0.95, that is, all coefficients of resti¬ 
tution for which the liquid-gas phase separation can be 
uniquely identified. 

Our simulations allow us to determine quantities that 
are not readily available experimentally, for example the 
pressure within the granular gas. We define the pressure 
in the homogeneous regions to be P = \{P XX + Pyy), 
where P xx and P yy are the horizontal components of the 
pressure tensor m■ To ensure we are in the steady state 
we first relax the system for ten seconds of simulated 
time. Thereafter pressure is averaged both spatially and 
over ten distinct initial configurations, for ten seconds 
each. To obtain the local pressure, the spatial average is 


additionally binned into boxes of size 5 d. 

By simulating small sample cells with horizontal di¬ 
mensions less than the liquid-gas interface width, phase 
separation can be suppressed. In this way the pres¬ 
sure can be calculated as a function of homogeneous 
quantities even under conditions for which a large sys¬ 
tem would phase separate. In our simulations, periodic 
boundary conditions are used in the horizontal directions. 
This method is a common practice for granular gases 
[2H EH (32H34], however, recent work questions whether 
the non-monotonic pressure-volume curves obtained rep¬ 
resent the equation of state for the system, or merely 
reflect finite size effects [35] [36]. Crucially we find that, 
for sufficiently small cells, the calculated pressure is not 
a function of the system size. 

Figure [2] shows the dependence of P on the dimension¬ 
less volume per particle, v = V/Nv p = in a small 
square-base cell of side L = 20 d for A=2.1d (solid line). 
As expected, the pressure exhibits a non-monotonic de¬ 
pendence on the volume, similar to what is observed in 
a molecular fluid. However, for the granular system, the 
pressure curve P(v ) is not an isotherm, and the physical 
origin of its non-monotonic shape is completely different. 
Unlike a molecular fluid, the granular gas has no attrac¬ 
tion between the particles; instead, the dilute phase is 
heated more effectively due to its intimate, resonant cou¬ 
pling to the vibrating walls, HUGH]. As a result, the 
non-monotonic behaviour in our system can be regarded 
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FIG. 3. Main panel: The solid lines show the pressures 
calculated in small cells for driving amplitudes in the range 
2.2d<A<2.8d. The filled circles show the corresponding bin- 
odal densities and the pressures calculated in long, thin sys¬ 
tems. The open circles show the pressure and densities pre¬ 
dicted using an equal-areas construction. Inset: Phase dia¬ 
gram for the liquid-gas-like phase separation. The filled and 
open circles show the binodal points determined by the long 
cell, and those predicted by P(v) and the equal-areas con¬ 
struction, respectively. The triangles show the spinodal points 
from the long cell and those predicted by the unstable region 
of P(v), respectively. 


as a crossover from a low temperature branch at high 
densities (left dashed curve in Fig. [2]) to a high tempera¬ 
ture branch at low densities (right dashed curve in Fig. [2]). 
The open symbols indicate the pressure calculated in cells 
of variable size, showing the convergence to the large-cell 
limit. We see that the region between the two extrema 
of P(v) is unstable against phase separation, resulting 
in a pressure corresponding to two-phase coexistence, P* 
(horizontal dashed line). 

Let us now turn to the equal-areas construction. We 
observe that the horizontal dashed line in Fig. [2] corre¬ 
sponding to P*, creates two approximately equal areas 
bounded above and below by the curve P(v) (hatched). 
To investigate this in more detail, the main panel of Fig. [3] 
shows P(v) calculated in the small cell for a range of driv¬ 
ing amplitudes. For each curve, the closed circles show 
the pressure determined using the large cell, the open 
circles indicate the pressure predicted by the equal-areas 
construction. It is quite striking that, for all the ampli¬ 
tudes investigated, the equal-areas construction predicts 
the pressure and the coexisting densities rather well. The 
inset in Fig. [3] shows the spinodal and binodal lines de¬ 
termined using the large cell (filled symbols), and the 
predictions made by using the small cell and equal-areas 
construction (open symbols). Again, we find remarkably 
good agreement given that we have no recourse to equi¬ 
librium thermodynamics, and as such there is no a-priori 



FIG. 4. The pressure-volume curve obtained from the small 
cell and an expanded view close to P*. The filled and open 
circles mark the averaged pressure and binodal volumes per 
particle for the liquid and gas phases respectively, obtained 
at different times. At each time, the average pressures in the 
two phases are equal to a good approximation. 


reason to expect the equal-areas construction to hold. 

In order to understand why the coexisting pressure, 
P* , in the phase-separated system is so close to the equal- 
areas pressure, Pea , obtained from the small-cell simula¬ 
tions, let us first consider a small cell with a fixed number 
of particles. If such a cell is expanded at some ambient 
pressure, from a volume corresponding to a mean liquid 
density <j> = (fn to that of the gas phase, <f> g , the work 
done per particle is readily obtained by integration of 
P{v), since P(v) does not depend on system size. If the 
expansion is carried out at Pea, then the total work done 
per particle vanishes. 

Consider now a large system with two coexisting 
phases. Mechanical equilibrium requires the pressure to 
be equal in both phases. Figure [4] shows the mean pres¬ 
sure, P, and volume per particle in each phase in the 
large cell, averaged over a few driving cycles. As the 
pressure fluctuates around P*, the corresponding densi¬ 
ties fluctuate so as to remain on the pressure curve P(v), 
confirming that P(v) does indeed act as the equation 
of state for our system. The pressure fluctuations are 
found to be tightly correlated with displacements of the 
interface between the two phases because the number of 
particles is fixed. These observations suggest the follow¬ 
ing simple model which is able to explain the origin of 
the equal-areas construction. 

Let the total volumes of the phases be V* = NiV t , where 
Ni and V\ are the number of particles and the specific 
volumes, respectively, for i € {l, g}. Consider a system of 
a fixed number of particles in a small box at a pressure P 
slightly away from P*, P = P*+V. The mechanical work 
done per particle to expand the box from specific volume 
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Vi to v* is, to leading order in V, w = ( v* — Vi)V. Here 
v* and v* are the corresponding liquid and gas specific 
volumes at P*. For any fluctuation, SNi + SN g = 0 and 
1 5Vi + SV g = 0. The volumes are related to the pressure 
through P(v) such that, to leading order, 5vt = —gtSV, 
where gi = —dv/dP \ v — v - is the compressibility. From 
these relations it follows that, for an additional small 
change in V , the change in particle number is SN g = 
[( giN t * + g g N*)/(v* — v*)]5V , where N* and N* are the 
number of liquid and gas particles at P*. Thus, the work 
done to move SN g particles across the interface is SWn = 
wSN g = (giN* + g g N*)V8V. For a finite fluctuation V 
the total work is 

Wtot = j^ 5 W n = i ( gi N? + g g N* g )V\ (1) 

which is quadratic in V. Consequently, any fluctuation 
that shifts the pressure away from P* while keeping the 
volumes per particle on P(v) requires mechanical work 
to be done. 

The energy needed to move the system away from P* 
is provided by momentary imbalances between the en¬ 
ergy injected by vibration and the energy dissipated in 
collisions. Because the system is dissipative, any resid¬ 
ual mechanical energy can always be released, moving 
the pressure back towards P*. This energy imbalance 
allows the bulk granular temperatures to fluctuate inde¬ 


pendently of the mechanical work. Although, in general, 
the equation of state for a granular system does depend 
on temperature, in a strongly vibrated system the tem¬ 
perature is fixed by the rate of energy injection and dis¬ 
sipation, Therefore, the fluctuations of the interface arise 
only from additional mechanical work done on one phase 
by the other. 

Although in a NESS the average dissipation balances 
with the injection, energy conservation is strongly vi¬ 
olated. Thus, it is quite unexpected that a global 
mechanical-energy minimization argument alone is suf¬ 
ficient to explain the behavior of the system. This is 
distinct from that of minimizing the local free energy, as 
applicable in thermodynamic equilibrium. In the main 
panel of Fig. [3] it is clear that the pressure calculated in 
the long cell deviates slightly but systematically from the 
equal-areas rule (only by less than 3%). It can be shown 
that this inequality follows from the higher order terms 
neglected in Eq. see supplemental material. 

Finally it is interesting to quote from Maxwell’s dis¬ 
cussion of the equal-areas rule in equilibrium systems: 
“Since the temperature has been constant throughout , no 
heat has been transformed into work” [25]. In our system 
the temperature is not constant throughout an expan¬ 
sion, yet to a good approximation an equal-areas rule still 
appears to be applicable based on the minimization of the 
residual mechanical work alone. It would be interesting 
to investigate whether this principle can be extended to 
a wider class of systems driven far from equilibrium. 
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